home *** CD-ROM | disk | FTP | other *** search
/ Aminet 6 / Aminet 6 - June 1995.iso / Aminet / gfx / 3d / irit50src.lha / irit5 / symb_lib / symbzero.c < prev   
Encoding:
C/C++ Source or Header  |  1995-02-14  |  15.4 KB  |  434 lines

  1. /******************************************************************************
  2. * SymbZero.c - computes the zeros and extremes of a given object.          *
  3. *******************************************************************************
  4. * Written by Gershon Elber, March 93.                          *
  5. ******************************************************************************/
  6.  
  7. #include "symb_loc.h"
  8.  
  9. #define SET_DEFAULT_EPSILON    1e-3
  10. #define SET_DEFAULT_EPSILON    1e-3
  11. #define ZERO_APX_EQ(x, y)        (ABS(x - y) < GlblSetEpsilon * 10)
  12. #define MID_RATIO            0.5123456789
  13.  
  14. static CagdPtStruct
  15.     *GlblPtList = NULL;
  16. static CagdRType CrvTMin, CrvTMax,
  17.     GlblSetEpsilon = SET_DEFAULT_EPSILON;
  18.  
  19. static void SymbScalarCrvConstSet(CagdCrvStruct *Crv, CagdRType ConstVal);
  20. static void SymbScalarCrvExtremSet(CagdCrvStruct *Crv);
  21. static void SymbInsertNewParam(CagdRType t);
  22.  
  23. /*****************************************************************************
  24. * DESCRIPTION:                                                               M
  25. * Computes the zero set of a given curve, in the given axis (1-3 for X-Z).   M
  26. *   Returned is a list of the zero set points holding the parameter values   M
  27. * att Pt[0] of each point.                             M
  28. *                                                                            *
  29. * PARAMETERS:                                                                M
  30. *   Crv:        To compute its zero set.                                     M
  31. *   Axis:       The axis of Crv to compute zero set for.                     M
  32. *   Epsilon:    Tolerance control.                                           M
  33. *                                                                            *
  34. * RETURN VALUE:                                                              M
  35. *   CagdPtStruct *:   List of parameter values form which Crv is zero in     M
  36. *                     axis Axis.                                             M
  37. *                                                                            *
  38. * KEYWORDS:                                                                  M
  39. *   SymbCrvZeroSet, zero set, symbolic computation                           M
  40. *****************************************************************************/
  41. CagdPtStruct *SymbCrvZeroSet(CagdCrvStruct *Crv, int Axis, CagdRType Epsilon)
  42. {
  43.     return SymbCrvConstSet(Crv, Axis, Epsilon, 0.0);
  44. }
  45.  
  46. /*****************************************************************************
  47. * DESCRIPTION:                                                               M
  48. * Computes the extremum set of a given curve, in given axis (1-3 for X-Z).   M
  49. *   Returned is a list of the extreme set points holding the parameter       M
  50. * values at Pt[0] of each point.                         M
  51. *   One could compute the derivative of the curve and find its zero set.     M
  52. *   However, for rational curves, this will double the degree and slow down  M
  53. * the computation considerably.                             M
  54. *                                                                            *
  55. * PARAMETERS:                                                                M
  56. *   Crv:        To compute its extremum set.                                 M
  57. *   Axis:       The axis of Crv to compute extremum set for.                 M
  58. *   Epsilon:    Tolerance control.                                           M
  59. *                                                                            *
  60. * RETURN VALUE:                                                              M
  61. *   CagdPtStruct *:   List of parameter values form which Crv has an         M
  62. *                     extremum value in axis Axis.                           M
  63. *                                                                            *
  64. * KEYWORDS:                                                                  M
  65. *   SymbCrvExtremSet, extremum set, symbolic computation                     M
  66. *****************************************************************************/
  67. CagdPtStruct *SymbCrvExtremSet(CagdCrvStruct *Crv, int Axis, CagdRType Epsilon)
  68. {
  69.     CagdCrvStruct *CrvW, *CrvX, *CrvY, *CrvZ, *TCrv,
  70.     *ScalarCrv = NULL;
  71.  
  72.     GlblSetEpsilon = Epsilon;
  73.  
  74.     SymbCrvSplitScalar(Crv, &CrvW, &CrvX, &CrvY, &CrvZ);
  75.  
  76.     switch (Axis) {
  77.     case 1:
  78.         if (CrvX)
  79.         ScalarCrv = CrvX;
  80.         else
  81.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  82.         break;
  83.     case 2:
  84.         if (CrvY)
  85.         ScalarCrv = CrvY;
  86.         else
  87.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  88.         break;
  89.     case 3:
  90.         if (CrvZ)
  91.         ScalarCrv = CrvZ;
  92.         else
  93.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  94.         break;
  95.     default:
  96.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  97.     }
  98.  
  99.     Crv = SymbCrvMergeScalar(CrvW, ScalarCrv, NULL, NULL);
  100.     if (CrvW)
  101.     CagdCrvFree(CrvW);
  102.     if (CrvX)
  103.     CagdCrvFree(CrvX);
  104.     if (CrvY)
  105.     CagdCrvFree(CrvY);
  106.     if (CrvZ)
  107.     CagdCrvFree(CrvZ);
  108.     
  109.     GlblPtList = NULL;
  110.     if (CAGD_IS_BEZIER_CRV(Crv)) {
  111.     /* We need to save domains. */
  112.     TCrv = CnvrtBezier2BsplineCrv(Crv);
  113.     CagdCrvFree(Crv);
  114.     Crv = TCrv;
  115.     }
  116.  
  117.     CagdCrvDomain(Crv, &CrvTMin, &CrvTMax);
  118.     SymbScalarCrvExtremSet(Crv);
  119.     CagdCrvFree(Crv);
  120.  
  121.     return GlblPtList;
  122. }
  123.  
  124. /*****************************************************************************
  125. * DESCRIPTION:                                                               M
  126. * Computes the constant set of a given curve, in the given axis          M
  127. * (1-3 for X-Z).                                 M
  128. *   Returned is a list of the constant set points holding the parameter         M
  129. * values at Pt[0] of each point.                         M
  130. *                                                                            *
  131. * PARAMETERS:                                                                M
  132. *   Crv:        To compute its constant set.                                 M
  133. *   Axis:       The axis of Crv to compute constant set for.                 M
  134. *   Epsilon:    Tolerance control.                                           M
  135. *   ConstVal:   The value at which to compute the constant set.             M
  136. *                                                                            *
  137. * RETURN VALUE:                                                              M
  138. *   CagdPtStruct *:   List of parameter values form which Crv has an         M
  139. *                     value of ConstVal in axis Axis.                        M
  140. *                                                                            *
  141. * KEYWORDS:                                                                  M
  142. *   SymbCrvConstSet, constant set, zero set, symbolic computation            M
  143. *****************************************************************************/
  144. CagdPtStruct *SymbCrvConstSet(CagdCrvStruct *Crv,
  145.                   int Axis,
  146.                   CagdRType Epsilon,
  147.                   CagdRType ConstVal)
  148. {
  149.     CagdCrvStruct *CrvW, *CrvX, *CrvY, *CrvZ, *TCrv,
  150.     *ScalarCrv = NULL;
  151.  
  152.     GlblSetEpsilon = Epsilon;
  153.  
  154.     SymbCrvSplitScalar(Crv, &CrvW, &CrvX, &CrvY, &CrvZ);
  155.  
  156.     switch (Axis) {
  157.     case 1:
  158.         if (CrvX)
  159.         ScalarCrv = CrvX;
  160.         else
  161.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  162.         break;
  163.     case 2:
  164.         if (CrvY)
  165.         ScalarCrv = CrvY;
  166.         else
  167.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  168.         break;
  169.     case 3:
  170.         if (CrvZ)
  171.         ScalarCrv = CrvZ;
  172.         else
  173.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  174.         break;
  175.     default:
  176.         SYMB_FATAL_ERROR(SYMB_ERR_OUT_OF_RANGE);
  177.     }
  178.  
  179.     Crv = SymbCrvMergeScalar(CrvW, ScalarCrv, NULL, NULL);
  180.     if (CrvW)
  181.     CagdCrvFree(CrvW);
  182.     if (CrvX)
  183.     CagdCrvFree(CrvX);
  184.     if (CrvY)
  185.     CagdCrvFree(CrvY);
  186.     if (CrvZ)
  187.     CagdCrvFree(CrvZ);
  188.     
  189.     GlblPtList = NULL;
  190.     if (CAGD_IS_BEZIER_CRV(Crv)) {
  191.     /* We need to save domains. */
  192.     TCrv = CnvrtBezier2BsplineCrv(Crv);
  193.     CagdCrvFree(Crv);
  194.     Crv = TCrv;
  195.     }
  196.  
  197.     CagdCrvDomain(Crv, &CrvTMin, &CrvTMax);
  198.     SymbScalarCrvConstSet(Crv, ConstVal);
  199.     CagdCrvFree(Crv);
  200.  
  201.     return GlblPtList;
  202. }
  203.  
  204. /*****************************************************************************
  205. * DESCRIPTION:                                                               *
  206. * Computes the zero set of a scalar curve. Curve might be rational.         *
  207. *   Assumes open end condition on the curve parametric domain.             *
  208. *   Zero set points are append to the GlblPtList point list.             *
  209. *                                                                            *
  210. * PARAMETERS:                                                                *
  211. *   Crv:        To compute its constant set.                                 *
  212. *   ConstVal:   Thep value at which to compute the constant set.         *
  213. *                                                                            *
  214. * RETURN VALUE:                                                              *
  215. *   void                                                                     *
  216. *****************************************************************************/
  217. static void SymbScalarCrvConstSet(CagdCrvStruct *Crv, CagdRType ConstVal)
  218. {
  219.     int i,
  220.     Len = Crv -> Length,
  221.         Sign = 0;
  222.     CagdRType
  223.     *ScalarPts = Crv -> Points[1],
  224.     *ScalarWPts = Crv -> Points[0];
  225.  
  226.     if (SymbCrvPosNegWeights(Crv)) {
  227.     i = 0;                  /* Force subdivision of the curve. */
  228.     }
  229.     else {
  230.     for (i = 0; i < Len; i++) {
  231.         CagdRType
  232.         V = (ScalarWPts ? ScalarPts[i] / ScalarWPts[i] :
  233.                                   ScalarPts[i]) - ConstVal;
  234.         int NewSign = SIGN(V);
  235.  
  236.         if (Sign * NewSign < 0)
  237.         break;
  238.         else if (NewSign)
  239.         Sign = NewSign;
  240.     }
  241.     }
  242.  
  243.     if (i < Len) {
  244.     CagdRType TMin, TMax, TMid;
  245.  
  246.     /* Control polygon is both positive and negative. */
  247.     CagdCrvDomain(Crv, &TMin, &TMax);
  248.     TMid = TMin * MID_RATIO + TMax * (1.0 - MID_RATIO);
  249.  
  250.     if (TMax - TMin < GlblSetEpsilon / 5.0) {    /* Small enough - stop. */
  251.         SymbInsertNewParam((TMax + TMin) / 2.0);
  252.     }
  253.     else {                          /* Needs to subdiv. */
  254.         CagdCrvStruct
  255.             *Crv1 = CagdCrvSubdivAtParam(Crv, TMid),
  256.         *Crv2 = Crv1 -> Pnext;
  257.  
  258.         SymbScalarCrvConstSet(Crv1, ConstVal);
  259.         SymbScalarCrvConstSet(Crv2, ConstVal);
  260.  
  261.         CagdCrvFree(Crv1);
  262.         CagdCrvFree(Crv2);
  263.     }
  264.     }
  265. }
  266.  
  267. /*****************************************************************************
  268. * DESCRIPTION:                                                               *
  269. * Computes the extrem set of a scalar curve. Curve might be rational.         *
  270. *   Assumes open end condition on the curve parametric domain.             *
  271. *   Extreme set points are append to the GlblPtList point list.             *
  272. *                                                                            *
  273. * PARAMETERS:                                                                *
  274. *   Crv:        To compute its extremum set.                                 *
  275. *                                                                            *
  276. * RETURN VALUE:                                                              *
  277. *   void                                                                     *
  278. *****************************************************************************/
  279. static void SymbScalarCrvExtremSet(CagdCrvStruct *Crv)
  280. {
  281.     int i,
  282.     Len = Crv -> Length,
  283.         Sign = 0;
  284.     CagdRType
  285.     *ScalarPts = Crv -> Points[1],
  286.     *ScalarWPts = Crv -> Points[0];
  287.  
  288.     if (SymbCrvPosNegWeights(Crv)) {
  289.     i = 0;                  /* Force subdivision of the curve. */
  290.     }
  291.     else {
  292.     CagdRType
  293.         PrevV = ScalarWPts ? ScalarPts[0] / ScalarWPts[0] : ScalarPts[0];
  294.  
  295.     for (i = 0; i < Len; i++) {
  296.         CagdRType
  297.         V = (ScalarWPts ? ScalarPts[i] / ScalarWPts[i] : ScalarPts[i]),
  298.         DV = V - PrevV;
  299.  
  300.         int NewSign = SIGN(DV);
  301.  
  302.         if (Sign * NewSign < 0)
  303.         break;
  304.         else if (NewSign)
  305.         Sign = NewSign;
  306.  
  307.         PrevV = V;
  308.     }
  309.     }
  310.  
  311.     if (i < Len) {
  312.     CagdRType TMin, TMax, TMid;
  313.  
  314.     /* Control polygon is both increasing and decreasing. */
  315.     CagdCrvDomain(Crv, &TMin, &TMax);
  316.     TMid = TMin * MID_RATIO + TMax * (1.0 - MID_RATIO);
  317.  
  318.     if (TMax - TMin < GlblSetEpsilon / 5.0) {    /* Small enough - stop. */
  319.         SymbInsertNewParam((TMax + TMin) / 2.0);
  320.     }
  321.     else {                          /* Needs to subdiv. */
  322.         CagdCrvStruct
  323.             *Crv1 = CagdCrvSubdivAtParam(Crv, TMid),
  324.         *Crv2 = Crv1 -> Pnext;
  325.  
  326.         SymbScalarCrvExtremSet(Crv1);
  327.         SymbScalarCrvExtremSet(Crv2);
  328.  
  329.         CagdCrvFree(Crv1);
  330.         CagdCrvFree(Crv2);
  331.     }
  332.     }
  333.     else {
  334.     CagdRType V1, V2, TMin, TMax;
  335.  
  336.     CagdCrvDomain(Crv, &TMin, &TMax);
  337.  
  338.     /* This segment is monotone. Test if end points are extremes. */
  339.     V1 = ScalarWPts ? ScalarPts[0] / ScalarWPts[0] : ScalarPts[0];
  340.     V2 = ScalarWPts ? ScalarPts[1] / ScalarWPts[1] : ScalarPts[1];
  341.     if (APX_EQ(V1, V2))
  342.         SymbInsertNewParam(TMin);
  343.  
  344.     V1 = ScalarWPts ? ScalarPts[Len - 2] / ScalarWPts[Len - 2]
  345.             : ScalarPts[Len - 2];
  346.     V2 = ScalarWPts ? ScalarPts[Len - 1] / ScalarWPts[Len - 1]
  347.             : ScalarPts[Len - 1];
  348.     if (APX_EQ(V1, V2))
  349.         SymbInsertNewParam(TMax);
  350.     }
  351. }
  352.  
  353. /*****************************************************************************
  354. * DESCRIPTION:                                                               M
  355. * Returns TRUE iff the Crv is not rational or rational with weights that are M
  356. * entirely positive or entirely negative.                     M
  357. *                                                                            *
  358. * PARAMETERS:                                                                M
  359. *   Crv:       To examine for same sign weights, if any.                     M
  360. *                                                                            *
  361. * RETURN VALUE:                                                              M
  362. *   CagdBType:   TRUE if no weights or all of same sign.                     M
  363. *                                                                            *
  364. * KEYWORDS:                                                                  M
  365. *   SymbCrvPosNegWeights, symbolic computation                               M
  366. *****************************************************************************/
  367. CagdBType SymbCrvPosNegWeights(CagdCrvStruct *Crv)
  368. {
  369.     int i;
  370.     CagdBType HasNeg, HasPos;
  371.     CagdRType
  372.     *Weights = Crv -> Points[0];
  373.  
  374.     if (Weights == NULL)
  375.     return FALSE;                   /* Curve is not rational. */
  376.  
  377.     for (HasNeg = HasPos = FALSE, i = Crv -> Length - 1; i >= 0; i--) {
  378.     HasNeg |= *Weights < 0.0;
  379.     HasPos |= *Weights++ > 0.0;
  380.     }
  381.  
  382.     return HasNeg && HasPos;
  383. }
  384.  
  385. /*****************************************************************************
  386. * DESCRIPTION:                                                               *
  387. * Insert a single t value into existing GlblPtList, provided no equal t      *
  388. * value exists already in the list. List is ordered incrementally.         *
  389. *                                                                            *
  390. * PARAMETERS:                                                                *
  391. *   t:         New value to insert to global GlblPtList list.             *
  392. *                                                                            *
  393. * RETURN VALUE:                                                              *
  394. *   void                                                                     *
  395. *****************************************************************************/
  396. static void SymbInsertNewParam(CagdRType t)
  397. {
  398.     CagdPtStruct *PtTmp, *PtLast, *Pt;
  399.  
  400.     if (APX_EQ(t, CrvTMin) || APX_EQ(t, CrvTMax))
  401.     return;
  402.  
  403.     Pt = CagdPtNew();
  404.     Pt -> Pt[0] = t;
  405.     Pt -> Pnext = NULL;
  406.  
  407.     if (GlblPtList) {
  408.     for (PtTmp = GlblPtList, PtLast = NULL;
  409.          PtTmp != NULL;
  410.          PtLast = PtTmp, PtTmp = PtTmp -> Pnext) {
  411.         if (ZERO_APX_EQ(PtTmp -> Pt[0], t)) {
  412.             IritFree(Pt);
  413.         return;
  414.         }
  415.         if (PtTmp -> Pt[0] > t)
  416.             break;
  417.     }
  418.     if (PtTmp) {
  419.         /* Insert the new point in the middle of the list. */
  420.         Pt -> Pnext = PtTmp;
  421.         if (PtLast)
  422.         PtLast -> Pnext = Pt;
  423.         else
  424.         GlblPtList = Pt;
  425.     }
  426.     else {
  427.         /* Insert the new point as the last point in the list. */
  428.         PtLast -> Pnext = Pt;
  429.     }
  430.     }
  431.     else
  432.         GlblPtList = Pt;
  433. }
  434.